parameters
if(is.null(input$cluster)) {
input$cluster = '0'
print("Warning: setting cluster to 0")
}
if(is.null(input$cluster5)) {
input$cluster5 = '0'
print("Warning: setting cluster5 to 0")
}
if(is.null(input$clusters)) {
input$clusters = '0'
print("Warning: setting clusters to 0")
}
## [1] "Warning: setting clusters to 0"
if(is.null(input$clusters1)) {
input$clusters1 = '0'
print("Warning: setting clusters1 to All")
}
if(is.null(input$clusters2)) {
input$clusters2 = '0'
print("Warning: setting clusters2 to 0")
}
if(is.null(input$clusters3)) {
input$clusters3 = '0'
print("Warning: setting clusters4 to 0")
}
if(is.null(input$clusters4)) {
input$clusters4 = 'All'
print("Warning: setting clusters4 to All")
}
dataTables = list()
dataTables$log2cpmOrg <- log2cpm
dataTables$tsne.data <- tsne.data
dataTables$tsne.dataOrg <- tsne.data
featuredata <-
featuredata[which(featuredata$Chromosome.Name %in% c(unlist(lapply(
seq(1, 22, 1), toString
)), c("X", "Y", "MT", "N"))),]
featuredata$Associated.Gene.Name <-
toupper(featuredata$Associated.Gene.Name)
featuredata <- featuredata[rownames(log2cpm), ]
dataTables$featuredataOrg <- featuredata
dataTables$positiveCells <- NULL
dataTables$positiveCellsAll <- NULL
# take only genes that are in all tables
rnames = rownames(featuredata)
rnames = rnames[rnames %in% rownames(log2cpm)]
rnames = rnames[rnames %in% rownames(gbm)]
dataTables$log2cpm <- log2cpm[rnames, ]
dataTables$gbm = gbm[rnames, ]
dataTables$gbm_log = gbm_log[rnames, ]
dataTables$featuredata <- featuredata[rnames, ]
geneC = colSums(log2cpm>0, na.rm = TRUE)
medianENSG = median(t(geneC))
umiC = colSums(log2cpm, na.rm = TRUE)
medianUMI = median(t(umiC))
minG = input$minGenes
gbm = as.matrix(exprs(dataTables$gbm))
goodCols = colSums(gbm) > minG
maxG = input$maxGenes
goodCols = goodCols & ( colSums(gbm) <= maxG)
geneNames = input$minExpGenes
ids = which(grepl(geneNames, dataTables$featuredata$Associated.Gene.Name))
goodCols = goodCols & (colSums(gbm[ids, ]) > 0)
useCells = goodCols
#genes to be explicitly removed
ipIDs = input$selectIds
if(nchar(ipIDs)>0) {
rmIds = !grepl(input$selectIds, dataTables$featuredata$Associated.Gene.Name)
}else
{
rmIds = rep(TRUE,nrow(dataTables$log2cpm))
}
# gene groups to be included
if(!is.null(input$geneListSelection)){
selectedgeneList = get_selected(input$geneListSelection)
if(length(selectedgeneList)>0){
selGenes = c()
for(sIdx in 1:length(selectedgeneList)){
print(sIdx)
att = attr(selectedgeneList[[sIdx]], "ancestry")
if(length(att)>0){
selGenes = c(selGenes , geneLists[[att]][[selectedgeneList[[sIdx]]]])
}
}
selGenes = unique(selGenes)
rmIds = rownames(dataTables$log2cpm) %in% selGenes & rmIds
}
}
# overall gene expression Min
minGene <- input$minGenesGS
if(!is.null(minGene)){
selGenes = rowSums(as.matrix(exprs(dataTables$gbm[,useCells]))) >=minGene
rmIds = rmIds & selGenes
}
if(DEBUG)cat(file=stderr(), "useGenes: done\n")
useGenes = rmIds
featureDataReact = dataTables$featuredata[useGenes, ]
featureData = dataTables$featuredata[useGenes, ]
if(DEBUG)cat(file=stderr(), "gbm\n")
gbm = dataTables$gbm[useGenes, useCells]
if(DEBUG)cat(file=stderr(), "gbm:DONE\n")
if(DEBUG)cat(file=stderr(), "gbm_log\n")
gbm_log = dataTables$gbm_log[useGenes, useCells]
if(DEBUG)cat(file=stderr(), "gbm_log:Done\n")
if(DEBUG)cat(file=stderr(), "scaterReads\n")
fd = featureDataReact
counts = as.matrix(exprs( gbm))
anno = pData(gbm)
anno$sample_id = anno$barcode
anno$fixed = "red"
# anno$individual= "NA1"
# anno$replicate = "r1"
# anno$well = "A01"
# anno$batch = "b1"
pheno_data <- new("AnnotatedDataFrame", anno)
rownames(pheno_data) <- pheno_data$sample_id
reads = as.matrix(counts)
rownames(reads) = make.unique(fd[rownames(reads),"Associated.Gene.Name"])
rownames(reads)[is.na(rownames(reads)) ] = "na"
reads = SingleCellExperiment(
assays = list(counts = reads),
colData = anno
)
ercc <- rownames(reads)[grepl("ERCC-", rownames(reads))]
mt = rownames(reads)[grepl("^MT",rownames(reads))]
reads <- scater::calculateQCMetrics(
reads
)
filter_by_expr_features <- (reads$total_features > 200)
reads$use <- (
# sufficient features (genes)
filter_by_expr_features
# sufficient molecules counted
#filter_by_total_counts &
# sufficient endogenous RNA
#filter_by_ERCC &
# remove cells with unusual number of reads in MT genes
#filter_by_MT
)
scaterReads = reads
if(DEBUG)cat(file=stderr(), "log2cpm\n")
if(DEBUG)cat(file=stderr(), "log2cpm:Done\n")
log2cpm = dataTables$log2cpm[useGenes, useCells]
if(DEBUG)cat(file=stderr(), "pca\n")
pca = run_pca(gbm_log)
## Variance explained by PCs: 0.02963652
if(DEBUG)cat(file=stderr(), "tsne\n")
set.seed(seed = seed)
tsne = run_tsne(pca, dims = 3, perplexity = 30, theta = 0.5)
## Using number of of pcs: 10
if(DEBUG)cat(file=stderr(), "tsne: done\n")
if(DEBUG)cat(file=stderr(), "kmClustering\n")
clustering=list()
kNr = 10
for(kNr in 2:10) {
set.seed(seed = seed)
km = run_kmeans_clustering(pca, k=kNr)
clustering[[paste0("kmeans_",kNr,"_clusters")]] = data.frame("Barcode" = rownames(data.frame(km$cluster)), "Cluster" = km$cluster)
}
## Using number of of pcs: 10
## Using number of of pcs: 10
## Using number of of pcs: 10
## Using number of of pcs: 10
## Using number of of pcs: 10
## Using number of of pcs: 10
## Using number of of pcs: 10
## Using number of of pcs: 10
## Using number of of pcs: 10
kmClustering = clustering
if(DEBUG)cat(file=stderr(), "kmClustering:done\n")
if(DEBUG)cat(file=stderr(), "tsne.data\n")
tsne.data = data.frame(tsne$Y)
colnames(tsne.data) = c("V1", "V2", "V3")
tsne.data$dbCluster = clustering$kmeans_10_clusters$Cluster-1
rownames(tsne.data) = clustering$kmeans_10_clusters$Barcode
if(DEBUG)cat(file=stderr(), "tsne.data: done\n")
set.seed(seed = seed)
prioritized_genes = prioritize_top_genes(gbm, tsne.data$dbCluster, "sseq",
logscale = FALSE,
min_mean=0.5,
p_cutoff=0.05,
order_by='pvalue')
## Computing differential expression parameters...
## Comparing Cluster 0 against the other clusters...
## Comparing Cluster 1 against the other clusters...
## Comparing Cluster 2 against the other clusters...
## Comparing Cluster 3 against the other clusters...
## Comparing Cluster 4 against the other clusters...
## Comparing Cluster 5 against the other clusters...
## Comparing Cluster 6 against the other clusters...
## Comparing Cluster 7 against the other clusters...
## Comparing Cluster 8 against the other clusters...
## Comparing Cluster 9 against the other clusters...
if(DEBUG)cat(file=stderr(), "dge\n")
if(!is.null(input$db1) & ! is.null(input$db2)){
subsetData <- subset(tsne.data, dbCluster == input$clusters1)
cells.1 <- rownames(brushedPoints(subsetData, input$db1))
cells.2 <- rownames(brushedPoints(subsetData, input$db2))
subsetExpression <- log2cpm[, union(cells.1, cells.2)]
genes.use <- rownames(subsetExpression)
data.1 = apply(subsetExpression[genes.use, cells.1], 1, expMean)
data.2 = apply(subsetExpression[genes.use, cells.2], 1, expMean)
total.diff = (data.1 - data.2)
genes.diff = names(which(abs(total.diff) > .2))
genes.use = ainb(genes.use, genes.diff)
toReturn <-
DiffExpTest(subsetExpression, cells.1, cells.2, genes.use = genes.use)
toReturn[, "avg_diff"] = total.diff[rownames(toReturn)]
toReturn$Associated.Gene.Name <-
featureData[rownames(toReturn), 'Associated.Gene.Name']
selectedDge <- toReturn
cat(stderr(), rownames(toReturn)[1:5])
dge <- toReturn
if(DEBUG)cat(file=stderr(), "dge:done11\n")
}else{
dge=NA
}
## 2 ENSG00000186891 ENSG00000116209 ENSG00000031698 ENSG00000265241 ENSG00000163191
# selectedDge <- reactiveValues()
p <-
plot_ly(
tsne.data,
x = ~ V1,
y = ~ V2,
z = ~ V3,
type = "scatter3d",
color = ~ dbCluster,
hoverinfo = "text",
text = paste('Cluster:', tsne.data$dbCluster),
mode = 'markers',
marker =
list(
line = list(width = 0),
size = rep(10, nrow(tsne.data)),
sizeref = 3
)
)
if(DEBUG)cat(file=stderr(), "output$tsne_main: done\n")
layout(p)
line1<-paste('No. of cells:', dim(log2cpm)[2],sep='\t')
line2<-paste('Median UMIs:', medianUMI,sep='\t')
line3<-paste('Median Genes:', medianENSG,sep='\t')
line5<-paste('No. of reads:', dim(log2cpm)[1],sep='\t')
HTML(
paste0("Summary statistics of this dataset:", '<br/>','<br/>',
line1, '<br/>',
line2, '<br/>',
line3, '<br/>',
line5
)
)
Summary statistics of this dataset:
No. of cells: 904
Median UMIs: 2231.12313456375
Median Genes: 1216
No. of reads: 14598
hist(colSums(as.matrix(exprs(gbm))), breaks = 50, main="histogram of number of UMI per cell")

barplot(pca$var_pcs, main="Variance captured by first PCs")

barplot(pca$var_pcs, main="Variance captured by first PCs")

gbmMat = as.matrix(exprs(dataTables$gbm))
fd = dataTables$featuredata
dt = fd[useGenes,c("Associated.Gene.Name", "Gene.Biotype", "Description")]
dt$rowSums = rowSums(gbmMat[useGenes,useCells])
dt$rowSamples = rowSums(gbmMat[useGenes,useCells]>0)
DT::datatable(dt)
## Warning in instance$preRenderHook(instance): It seems your data is too
## big for client-side DataTables. You may consider server-side processing:
## https://rstudio.github.io/DT/server.html
gbmMat = as.matrix(exprs(dataTables$gbm))
fd = dataTables$featuredata
dt = fd[useGenes,c("Associated.Gene.Name", "Gene.Biotype", "Description")]
dt$rowSums = rowSums(gbmMat[useGenes,useCells])
dt$rowSamples = rowSums(gbmMat[useGenes,useCells]>0)
DT::datatable(dt)
## Warning in instance$preRenderHook(instance): It seems your data is too
## big for client-side DataTables. You may consider server-side processing:
## https://rstudio.github.io/DT/server.html
if(DEBUG)cat(file=stderr(), "output$tsne_plt\n")
geneid <- rownames(featureData[which(featureData$Associated.Gene.Name ==
toupper(input$gene_id)), ])[1]
expression <- log2cpm[geneid, ]
tsneData <- cbind(tsne.data, t(expression))
names(tsneData)[names(tsneData) == geneid] <- 'values'
p <-
plot_ly(
tsneData,
x = ~ V1,
y = ~ V2,
z = ~ V3,
type = "scatter3d",
hoverinfo = "text",
text = paste('Cluster:', tsne.data$dbCluster),
mode = 'markers',
marker = list(
size = 2,
line = list(width = 0),
color = ~ values,
colors = 'Greens'
)
)
layout(p, title = toupper(input$gene_id))
if(DEBUG)cat(file=stderr(), "output$clusterPlot\n")
geneid <- rownames(featureData[which(featureData$Associated.Gene.Name ==
toupper(input$gene_id)), ])[1]
expression <- log2cpm[geneid, ]
tsneData <- cbind(tsne.data, t(expression))
names(tsneData)[names(tsneData) == geneid] <- 'values'
if(is.null(input$cluster)) input$cluster = 0
subsetData <- subset(tsneData, dbCluster == input$cluster)
## Warning in dbCluster == input$cluster: longer object length is not a
## multiple of shorter object length
p1 <-
ggplot(subsetData,
aes_string(x = input$dimension_x, y = input$dimension_y)) +
geom_point(aes_string(size = 2, color = 'values')) +
geom_point(shape = 1,
size = 4,
colour = "black") +
theme_bw() +
theme(
axis.text.x = element_text(
angle = 90,
size = 12,
vjust = 0.5
),
axis.text.y = element_text(size = 12),
strip.text.x = element_text(size = 16),
strip.text.y = element_text(size = 14),
axis.title.x = element_text(face = "bold", size = 16),
axis.title.y = element_text(face = "bold", size = 16),
legend.position = "none"
) +
ggtitle(paste(toupper(input$gene_id), input$cluster, sep = '-Cluster')) +
scale_colour_gradient2(low = 'grey50', high = "red")
p1

if(DEBUG)cat(file=stderr(), "output$gene_vio_plot\n")
geneid <- rownames(featureData[which(featureData$Associated.Gene.Name ==
toupper(input$gene_id)), ])[1]
expression <- log2cpm[geneid, ]
validate(need(is.na(sum(expression)) != TRUE, ''))
tsneData <- cbind(tsne.data, t(expression))
names(tsneData)[names(tsneData) == geneid] <- 'values'
#tsne.data<-subset(tsne.data,dbCluster!=0)
p1 <-
ggplot(tsneData, aes(factor(dbCluster), values, fill = factor(dbCluster))) +
geom_violin(scale = "width") +
stat_summary(
fun.y = median,
geom = "point",
size = 5,
color = 'black'
) +
stat_summary(fun.data = n_fun, geom = "text") +
theme_bw() +
theme(
axis.text.x = element_text(
angle = 90,
size = 12,
vjust = 0.5
),
axis.text.y = element_text(size = 12),
strip.text.x = element_text(size = 16),
strip.text.y = element_text(size = 14),
axis.title.x = element_text(face = "bold", size = 16),
axis.title.y = element_text(face = "bold", size = 16),
legend.position = "none"
) +
xlab('Cluster') +
ylab('Expression') +
ggtitle(toupper(input$gene_id))
if(DEBUG)cat(file=stderr(), "output$gene_vio_plot:done\n")
p1

if(DEBUG)cat(file=stderr(), "output$panelPlot\n")
genesin <- input$panelplotids
genesin <- toupper(genesin)
genesin <- gsub(" ", "", genesin, fixed = TRUE)
genesin <- strsplit(genesin, ',')
genesin<-genesin[[1]]
par(mfrow=c(ceiling(length(genesin)/4),4), mai = c(0, 0., 0., 0.))
rbPal <- colorRampPalette(c('#f0f0f0','red'))
if (input$clusters4 == 'All')
{
for (i in 1:length(genesin)){
Col <- rbPal(10)[
as.numeric(
cut(
as.numeric(
log2cpm[
rownames(featureData[which(featureData$Associated.Gene.Name==genesin[i]),])
,]
),breaks = 10))]
plot(tsne.data[,input$dimension_x4],tsne.data[,input$dimension_y4],col=Col,pch=16,axes = FALSE,frame.plot = TRUE, ann=FALSE)
title(genesin[i],line=-1.2,adj = 0.05,cex.main=2)
}
}else{
for (i in 1:length(genesin)){
subsetTSNE <- subset(tsne.data, dbCluster == input$clusters4)
Col <- rbPal(10)[
as.numeric(
cut(
as.numeric(
log2cpm[
rownames(featureData[which(featureData$Associated.Gene.Name==genesin[i]),])
,]
),breaks = 10))]
names(Col)<-rownames(tsne.data)
plotCol<-Col[rownames(subsetTSNE)]
plot(subsetTSNE[,input$dimension_x4],subsetTSNE[,input$dimension_y4],col=plotCol,pch=16,axes = FALSE,frame.plot = TRUE, ann=FALSE)
title(genesin[i],line=-1.2,adj = 0.05,cex.main=2)
}
}

scater::plotQC(scaterReads, type = "highest-expression", col_by_variable="fixed")

genesin <- input$heatmap_geneids
genesin <- toupper(genesin)
genesin <- gsub(" ", "", genesin, fixed = TRUE)
genesin <- strsplit(genesin, ',')
map <- rownames(featureData[which(featureData$Associated.Gene.Name %in% genesin[[1]]), ])
cat(file = stderr(), length(map))
expression <- log2cpm[map, ]
if(is.na(sum(expression)) == TRUE){
print('Gene symbol incorrect or genes not expressed')
}else{
tsne.data <- tsne.data[order(tsne.data$dbCluster), ]
expression <- expression[, rownames(tsne.data)]
expression <- expression[complete.cases(expression), ]
annotation <- data.frame(factor(tsne.data$dbCluster))
rownames(annotation) <- colnames(expression)
colnames(annotation) <- c('Cluster')
h <-
pheatmap(
as.matrix(expression),
cluster_rows = TRUE,
cluster_cols = FALSE,
scale = 'row',
fontsize_row = 10,
labels_col = colnames(expression),
labels_row = featureData[rownames(expression), 'Associated.Gene.Name'],
show_rownames = TRUE,
annotation_col = annotation,
show_colnames = FALSE,
annotation_legend = TRUE,
breaks = seq(-6, 6, by = .12),
colorRampPalette(rev(brewer.pal(
n = 6, name =
"RdBu"
)))(100)
)
print(h)
}

## $tree_row
##
## Call:
## hclust(d = d, method = method)
##
## Cluster method : complete
## Distance : euclidean
## Number of objects: 4
##
##
## $tree_col
## [1] NA
##
## $kmeans
## [1] NA
##
## $gtable
## TableGrob (5 x 6) "layout": 7 grobs
## z cells name grob
## 1 1 (4-4,1-1) row_tree polyline[GRID.polyline.16057]
## 2 2 (4-4,3-3) matrix gTree[GRID.gTree.16059]
## 3 3 (4-4,4-4) row_names text[GRID.text.16060]
## 4 4 (3-3,3-3) col_annotation rect[GRID.rect.16061]
## 5 5 (3-3,4-4) col_annotation_names text[GRID.text.16062]
## 6 6 (3-5,6-6) annotation_legend gTree[GRID.gTree.16067]
## 7 7 (3-5,5-5) legend gTree[GRID.gTree.16070]
if(DEBUG)cat(file=stderr(), "output$clusterPlot2\n")
geneid <- rownames(featureData[which(featureData$Associated.Gene.Name ==
toupper(input$gene_id_sch)), ])[1]
expression <- log2cpm[geneid, ]
validate(need(is.na(sum(expression)) != TRUE, ''))
tsneData <- cbind(tsne.data, t(expression))
names(tsneData)[names(tsneData) == geneid] <- 'values'
# if(DEBUG)cat(file=stderr(), paste("output$dge_plot1:---",input$clusters2,"---\n"))
if(!is.null(input$dimension_x2) & !is.null(input$dimension_y2)){
subsetData <- subset(tsneData, dbCluster == input$clusters2)
p1 <-
ggplot(subsetData,
aes_string(x = input$dimension_x2, y = input$dimension_y2)) +
geom_point(aes_string(size = 2, color = 'values')) +
geom_point(shape = 1,
size = 4,
colour = "black") +
theme_bw() +
theme(
axis.text.x = element_text(
angle = 90,
size = 12,
vjust = 0.5
),
axis.text.y = element_text(size = 10),
strip.text.x = element_text(size = 16),
strip.text.y = element_text(size = 14),
axis.title.x = element_text(face = "bold", size = 16),
axis.title.y = element_text(face = "bold", size = 16),
legend.position = "none"
) +
ggtitle(paste(toupper(input$gene_id_sch), input$clusters2, sep =
'-Cluster')) +
scale_colour_gradient2(low = 'grey50', high = "red")
print(p1)
}

if(DEBUG)cat(file=stderr(), "output$selectedHeatmap\n")
genesin <- input$heatmap_geneids2
genesin <- toupper(genesin)
genesin <- strsplit(genesin, ',')
if(!is.null(input$scb1)){
subsetData <-
subset(tsne.data, dbCluster == input$clusters2)
cells.1 <- rownames(brushedPoints(df=subsetData,
brush=input$scb1,
xvar = input$scb1$mapping$x,
yvar = input$scb1$mapping$y,
panelvar1 =
))
map <- rownames(featureData[which(featureData$Associated.Gene.Name %in% genesin[[1]]), ])
expression <- log2cpm[map, cells.1]
expression <- expression[complete.cases(expression), ]
cat(file = stderr(), rownames(expression))
mColor <- max(expression)
if( is.na(sum(expression)) == TRUE){
print('Gene symbol incorrect or genes not expressed')
}else{
h <-
pheatmap(
as.matrix(expression),
cluster_rows = TRUE,
cluster_cols = TRUE,
scale = 'row',
fontsize_row = 10,
labels_col = colnames(expression),
labels_row = featureData[rownames(expression), 'Associated.Gene.Name'],
show_rownames = TRUE,
show_colnames = FALSE,
breaks = seq(-6, 6, by = .12),
colorRampPalette(rev(brewer.pal(
n = 6, name =
"RdBu"
)))(100)
)
print(h)
}
}

## $tree_row
##
## Call:
## hclust(d = d, method = method)
##
## Cluster method : complete
## Distance : euclidean
## Number of objects: 2
##
##
## $tree_col
##
## Call:
## hclust(d = d, method = method)
##
## Cluster method : complete
## Distance : euclidean
## Number of objects: 44
##
##
## $kmeans
## [1] NA
##
## $gtable
## TableGrob (5 x 6) "layout": 5 grobs
## z cells name grob
## 1 1 (2-2,3-3) col_tree polyline[GRID.polyline.16121]
## 2 2 (4-4,1-1) row_tree polyline[GRID.polyline.16122]
## 3 3 (4-4,3-3) matrix gTree[GRID.gTree.16124]
## 4 4 (4-4,4-4) row_names text[GRID.text.16125]
## 5 5 (3-5,5-5) legend gTree[GRID.gTree.16128]
if(DEBUG)cat(file=stderr(), "output$plotCoExpression\n")
genesin <- input$mclustids
genesin <- toupper(genesin)
genesin <- strsplit(genesin, ',')
subsetData <-
subset(tsne.data, dbCluster == input$clusters3)
cells.1 <- rownames(subsetData)
map <- rownames(featureData[which(featureData$Associated.Gene.Name %in% genesin[[1]]), ])
expression <- log2cpm[map, ]
if(is.na(sum(expression)) == TRUE){
print ('Gene symbol incorrect or genes not expressed')
}else{
bin <- expression
bin[] <- 0
for (i in 1:nrow(expression))
{
x <- Mclust(expression[i, ], G = 2)
bin[i, ] <- x$classification
}
bin <- bin - 1
allexprs <- apply(bin, 2, sum)
plotexprs <- allexprs
plotexprs[] <- 0
plotexprs[allexprs >= length(rownames(bin))] <- 1
positiveCells <- allexprs >= length(rownames(bin))
positiveCellsAll <- plotexprs
mergeExprs <- plotexprs[rownames(subsetData)]
subsetData$CoExpression <- mergeExprs
p1 <-
ggplot(subsetData,
aes_string(x = input$dimension_x3, y = input$dimension_y3)) +
geom_point(aes_string(size = 2, color = 'CoExpression')) +
geom_point(shape = 1,
size = 4,
colour = "black") +
theme_bw() +
theme(
axis.text.x = element_text(
angle = 90,
size = 12,
vjust = 0.5
),
axis.text.y = element_text(size = 12),
strip.text.x = element_text(size = 16),
strip.text.y = element_text(size = 14),
axis.title.x = element_text(face = "bold", size = 16),
axis.title.y = element_text(face = "bold", size = 16),
legend.position = "none"
) +
#ggtitle(paste(toupper(input$gene_id),input$cluster,sep='-Cluster'))+
scale_colour_gradient2(low = 'grey50', high = "red")
print(p1)
}
## fitting ...
##
|
| | 0%
|
|====================== | 33%
|
|=========================================== | 67%
|
|=================================================================| 100%

if(DEBUG)cat(file=stderr(), "output$onOffTable\n")
merge <- tsne.data
# if(DEBUG)cat(file=stderr(), paste("positiveCellsAll:---",positiveCellsAll,"---\n"))
merge$CoExpression <- positiveCellsAll
df <-
as.data.frame(table(merge[, c('dbCluster', 'CoExpression')]))
dfOut <- cast(df, dbCluster ~ CoExpression)
## Using Freq as value column. Use the value argument to cast to override this choice
colnames(dfOut) <- c("Cluster", 'OFF', 'ON')
rownames(dfOut) <- dfOut$Cluster
dfOut['Sum', ] <- c('', sum(dfOut$OFF), sum(dfOut$ON))
## Warning in `[<-.factor`(`*tmp*`, iseq, value = ""): invalid factor level,
## NA generated
DT::datatable(dfOut)
if(!is.null(input$dimension_x1) & !is.null(input$dimension_y1)){
subsetData <- subset(tsne.data, dbCluster == input$clusters1)
p1 <-
ggplot(subsetData,
aes_string(x = input$dimension_x1, y = input$dimension_y1)) +
geom_point() +
geom_point(shape = 1,
size = 4,
color = "black") +
theme_bw() +
theme(
axis.text.x = element_text(
angle = 90,
size = 12,
vjust = 0.5
),
axis.text.y = element_text(size = 12),
strip.text.x = element_text(size = 16),
strip.text.y = element_text(size = 14),
axis.title.x = element_text(face = "bold", size = 16),
axis.title.y = element_text(face = "bold", size = 16),
legend.position = "none"
) +
ggtitle(input$clusters1)
print(p1)
}

if(DEBUG)cat(file=stderr(), "output$dge_plot2\n")
if(!is.null(input$dimension_x1) & !is.null(input$dimension_y1)){
subsetData <- subset(tsne.data, dbCluster == input$clusters1)
p1 <-
ggplot(subsetData,
aes_string(x = input$dimension_x1, y = input$dimension_y1)) +
geom_point() +
geom_point(shape = 1,
size = 4,
color = "black") +
theme_bw() +
theme(
axis.text.x = element_text(
angle = 90,
size = 12,
vjust = 0.5
),
axis.text.y = element_text(size = 12),
strip.text.x = element_text(size = 16),
strip.text.y = element_text(size = 14),
axis.title.x = element_text(face = "bold", size = 16),
axis.title.y = element_text(face = "bold", size = 16),
legend.position = "none"
) +
ggtitle(input$clusters1)
print(p1)
}

if(DEBUG)cat(file=stderr(), "output$dge\n")
if(!is.na(dge)){
top.genes <- dge
top.genes$Associated.Gene.Name <-
featureData[rownames(top.genes), 'Associated.Gene.Name']
if (dim(top.genes)[1] > 1) {
DT::datatable(top.genes,
options = list(
orderClasses = TRUE,
lengthMenu = c(10, 30, 50),
pageLength = 10
))
}
}
## Warning in if (!is.na(dge)) {: the condition has length > 1 and only the
## first element will be used
if(DEBUG)cat(file=stderr(), "output$crHeat_plot1\n")
example_K <- 10
example_Cols <- rev(brewer.pal(10,"Set3")) # customize plotting colors
cells_to_plot <- order_cell_by_clusters(gbm, tsne.data$dbCluster)
example_col = example_Cols[1:example_K]
p = gbm_pheatmap(gbm=log_gene_bc_matrix(gbm),
genes_to_plot=prioritized_genes,
cells_to_plot=cells_to_plot,
n_genes=10,
colour=example_col,
limits = c(-3, 3))

if(DEBUG)cat(file=stderr(), "output$crHeat_plot1:done\n")
p
## $tree_row
## [1] NA
##
## $tree_col
## [1] NA
##
## $kmeans
## [1] NA
##
## $gtable
## TableGrob (5 x 6) "layout": 6 grobs
## z cells name grob
## 1 1 (4-4,3-3) matrix gTree[GRID.gTree.16274]
## 2 2 (4-4,4-4) row_names text[GRID.text.16275]
## 3 3 (3-3,3-3) col_annotation rect[GRID.rect.16276]
## 4 4 (4-4,2-2) row_annotation rect[GRID.rect.16277]
## 5 5 (3-5,6-6) annotation_legend gTree[GRID.gTree.16285]
## 6 6 (3-5,5-5) legend gTree[GRID.gTree.16288]
if(DEBUG)cat(file=stderr(), "output$crPrioGenes\n")
dt = prioritized_genes[[input$cluster5]]
DT::datatable(dt)
## Warning in instance$preRenderHook(instance): It seems your data is too
## big for client-side DataTables. You may consider server-side processing:
## https://rstudio.github.io/DT/server.html
print(ls())
## [1] "Col" "DEBUG"
## [3] "allexprs" "anno"
## [5] "annotation" "bin"
## [7] "cells.1" "cells.2"
## [9] "cells_to_plot" "clustering"
## [11] "counts" "data.1"
## [13] "data.2" "dataTables"
## [15] "df" "dfOut"
## [17] "dge" "dt"
## [19] "ercc" "example_Cols"
## [21] "example_K" "example_col"
## [23] "expression" "fd"
## [25] "featureData" "featureDataReact"
## [27] "featuredata" "filter_by_expr_features"
## [29] "gbm" "gbmMat"
## [31] "gbm_log" "geneC"
## [33] "geneNames" "geneid"
## [35] "genes.diff" "genes.use"
## [37] "genesin" "goodCols"
## [39] "h" "i"
## [41] "ids" "input"
## [43] "ipIDs" "kNr"
## [45] "km" "kmClustering"
## [47] "line1" "line2"
## [49] "line3" "line5"
## [51] "log2cpm" "mColor"
## [53] "map" "maxG"
## [55] "medianENSG" "medianUMI"
## [57] "merge" "mergeExprs"
## [59] "minG" "minGene"
## [61] "mt" "p"
## [63] "p1" "params"
## [65] "pca" "pheno_data"
## [67] "plotexprs" "positiveCells"
## [69] "positiveCellsAll" "prioritized_genes"
## [71] "rbPal" "reads"
## [73] "rmIds" "rnames"
## [75] "scaterReads" "seed"
## [77] "selGenes" "selectedDge"
## [79] "selectedgeneList" "subsetData"
## [81] "subsetExpression" "toReturn"
## [83] "top.genes" "total.diff"
## [85] "tsne" "tsne.data"
## [87] "tsneData" "umiC"
## [89] "useCells" "useGenes"
## [91] "x"
sessionInfo()
## R version 3.4.3 (2017-11-30)
## Platform: x86_64-apple-darwin15.6.0 (64-bit)
## Running under: macOS Sierra 10.12.6
##
## Matrix products: default
## BLAS: /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/3.4/Resources/lib/libRlapack.dylib
##
## locale:
## [1] C
##
## attached base packages:
## [1] stats4 parallel stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] bindrcpp_0.2 scater_1.6.3
## [3] SingleCellExperiment_1.0.0 SummarizedExperiment_1.8.1
## [5] DelayedArray_0.4.1 matrixStats_0.53.1
## [7] GenomicRanges_1.30.3 GenomeInfoDb_1.14.0
## [9] IRanges_2.12.0 S4Vectors_0.16.0
## [11] shinyWidgets_0.4.1 kableExtra_0.7.0
## [13] knitr_1.20 SCORPIUS_1.0
## [15] cellrangerRkit_2.0.0 Rmisc_1.5
## [17] plyr_1.8.4 lattice_0.20-35
## [19] bit64_0.9-7 bit_1.1-12
## [21] Biobase_2.38.0 BiocGenerics_0.24.0
## [23] Matrix_1.2-12 reshape_0.8.7
## [25] mclust_5.4 RColorBrewer_1.1-2
## [27] sm_2.2-5.4 shinyTree_0.2.2
## [29] threejs_0.3.1 igraph_1.1.2
## [31] pheatmap_1.0.8 edgeR_3.20.8
## [33] limma_3.34.9 DT_0.4
## [35] shinythemes_1.1.1 plotly_4.7.1
## [37] ggplot2_2.2.1 shinydashboard_0.6.1
## [39] shiny_1.0.5
##
## loaded via a namespace (and not attached):
## [1] Rtsne_0.13 ggbeeswarm_0.6.0 colorspace_1.3-2
## [4] rjson_0.2.15 rprojroot_1.3-2 XVector_0.18.0
## [7] base64enc_0.1-3 AnnotationDbi_1.40.0 xml2_1.2.0
## [10] ranger_0.9.0 codetools_0.2-15 splines_3.4.3
## [13] tximport_1.6.0 jsonlite_1.5 readr_1.1.1
## [16] compiler_3.4.3 httr_1.3.1 backports_1.1.2
## [19] assertthat_0.2.0 lazyeval_0.2.1 htmltools_0.3.6
## [22] prettyunits_1.0.2 tools_3.4.3 gtable_0.2.0
## [25] glue_1.2.0 GenomeInfoDbData_1.0.0 reshape2_1.4.3
## [28] dplyr_0.7.4 Rcpp_0.12.15 iterators_1.0.9
## [31] crosstalk_1.0.0 stringr_1.3.0 rvest_0.3.2
## [34] mime_0.5 irlba_2.3.2 XML_3.98-1.10
## [37] princurve_1.1-12 zlibbioc_1.24.0 MASS_7.3-49
## [40] scales_0.5.0 TSP_1.1-5 hms_0.4.1
## [43] rhdf5_2.22.0 yaml_2.1.16 memoise_1.1.0
## [46] pbapply_1.3-4 gridExtra_2.3 biomaRt_2.34.2
## [49] stringi_1.1.6 RSQLite_2.0 foreach_1.4.4
## [52] rlang_0.2.0 pkgconfig_2.0.1 bitops_1.0-6
## [55] evaluate_0.10.1 purrr_0.2.4 bindr_0.1
## [58] labeling_0.3 htmlwidgets_1.0 magrittr_1.5
## [61] R6_2.2.2 DBI_0.7 pillar_1.2.0
## [64] fitdistrplus_1.0-9 survival_2.41-3 RCurl_1.95-4.10
## [67] tibble_1.4.2 rmarkdown_1.8 viridis_0.5.0
## [70] progress_1.1.2 locfit_1.5-9.1 grid_3.4.3
## [73] data.table_1.10.4-3 blob_1.1.0 digest_0.6.15
## [76] xtable_1.8-2 tidyr_0.8.0 httpuv_1.3.5
## [79] munsell_0.4.3 beeswarm_0.2.3 viridisLite_0.3.0
## [82] vipor_0.4.5